We can simulate strings, membranes, and plates in different ways. For the linear case, it is practical to use a fast linear solver, such as solve_sinusoidal, or a an iir filter (iir_filter_parallel). The latter being slower than the former.
For most cases the steps are similar:
Define the parameters of the simulation
Get the eigenpairs
Discretize
Simulate
Linear models
For the linear models we want to solve the following equation:
n_modes =50n_steps =44100sample_rate =44100dt =1.0/ sample_rateexcitation_position =0.2readout_position =0.5initial_deflection =0.03n_gridpoints =101# number of gridpoints for evaluating the eigenfunctionsstring_params = StringParameters()
Get \(\gamma_{\mu}\) and \(\omega_{\mu}\) and integrate in time from the initial conditions defined above. This should be very fast, even for a large number of modes.
Transform the modal solution back to the physical space using the precomputed eigenfunctions, or evaluate at a single point.
mu = np.arange(1, n_modes +1) # mode indicesreadout_weights = evaluate_string_eigenfunctions( mu, readout_position, string_params,)# at a single pointu_readout = readout_weights @ modal_sol# at all pointssol = inverse_STL(K, modal_sol, string_params.length)
Single point:
Code
display(Audio(u_readout, rate=sample_rate))fig, ax = plt.subplots(1, 1, figsize=(6, 2))ax.plot(u_readout)ax.set_xlabel("Sample")ax.set_ylabel("Deflection (m)")ax.set_title("Deflection of the string at a single point")ax.set_xlim(-2, sample_rate //10)ax.grid(True)
All points:
Code
n_frames =2000frame_stride =15positions = np.linspace(0, string_params.length, sol.shape[0])time_samples = sol[:, :n_frames]fig, ax = plt.subplots(figsize=(6, 3))(line,) = ax.plot(positions, time_samples[:, 0])ax.set_ylim(np.min(time_samples) *1.1, np.max(time_samples) *1.1)ax.set_xlabel("Position (m)")ax.set_ylabel("Deflection (m)")ax.set_title("String Vibration")ax.grid(True)fig.tight_layout()def update(i): frame = i * frame_stride line.set_ydata(time_samples[:, frame])return (line,)ani = animation.FuncAnimation( fig, update, frames=n_frames // frame_stride, interval=25, blit=True,)# Only display the animation, not the still figureplt.close(fig)HTML(ani.to_jshtml())
Plate
We do a very similar thing for plates. First define the parameters